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Abstract 


This thesis is concerned with the development of analytical methods 
for the simulation of the complete ballistocardiogram (BCG). 

The first chapter of this thesis contains a derivation of the 
basic equation of BCG, definition of components of BCG which can be treated 
individually, and a critical note on BCG simulation as used by other authors. 

The second chapter presents a new method of simulation of the 
contribution of arterial circulation to the BCG. This method based on the 
assumption of linearity and using the frequency domain, is formulated as a 
digital computer program. Results obtained from this program can Aceiee 
and economically be used for BCG simulation. The method is accurate and can 
be used for the investigation of various linear models of arterial circulation 
or as a diagnostic tool for the design of analogs. Two of the results 
obtained in this chapter are : Modeling of the system as having distributed 
*pacameters, A new representation of small vessels for the purpose of BCG 
simulation, more realistic than models used by analogs and easy to handle 
by digital computers. 

The third chapter contains an analysis of factors causing the motion 
of the heart. It is more complete and realistic than known methods and it 
is completely analytical (it requires knowledge of mechanical properties of 
the tissue surrounding the heart and some hemodynamic variables). Some 


of the results obtained in this chapter agree well with reality, one factor, 
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however, does not seem to be represented adequately. It is suggested 
that better knowledge of some hemodynamic variables and a more accurate 
geometrical model account for discrepancies between these results and 
reality. 

The conclusion of the research presented here is that the 
description of components of the BCG is now complete and their simulation 


possible. 
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Introduction 


Ballistocardiography is a recording technique based on the 
fact that during the heart cycle, the distribution of mass in the human 
body is changing. At present it is generally accepted that the best way 
in which to make use of this phenomenon is to measure the displacement 
velocity or neeel eration of a bed supporting the subject in such a way 
that the subject is practically freely floating with the bed without much 
restraining mechanical coupling to the support of the bed. This technique 
is called ultralow frequency ballistocardiography. 

Ballistocardiography (BCG) is an old noninvasive method. It 
gives a quantitative picture of the pumping function of the heart and 
the resulting circulation - as opposed to electrocardiography, (ECG), 
which measures the electrical phenomena associated with the heart cycle. 
This is a good reason why BCG should be a widely accepted clinical 
method - to the same extent as ECG. In reality, however, relatively few 
hospitals use BCG, while ECG is generally accepted. One of the main 
reasons why this is so is that the BCG record (BCG) is not so well 
understood as the ECG. This is partly because the shape of BCG is more 
complex than ECG and its variability is larger, and partly because BCG 
represents the integrated effect of heart function on the whole 
circulatory system and individual peaks of the record are, therefore, not 
immediately associable with definable single causes. This is the reason 


why a relatively large number of researchers have in the last ten years 
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been attempting BCG simulation by various methods. Their work has 
already given many useful results but is far from being complete. The 
best known results have been obtained by Dr. Noordergraaf and his 
students from the Department of Biomedical Engineering at the University 
of Pennsylvania. Their approach is simulation of the human circulation 
by an electrical analog adapted to BCG simulation. It is described in (1, 
2, 3). (In the rest of this chapter analog A always means this particular 
analog). Analog A, however, cannot simulate all components of the BCG 
and gives only what its designers agree is the essential part of BCG. 
This fact has lead the author of this thesis to an investigation of an 
important component of the BCG not simulated by A - the contribution of 
the heart itself. This investigation is described in the third chapter 
of this thesis. Before the study of this component was undertaken, the 
basic equations of BCG simulation were reformulated and various components 
of BCG defined. This forms the first, chapter of this thesis.. The second 
chapter examines some of these components - those simulated by A. A new 
method of simulation is suggested - using a digital rather than analog 
computer. Since this approach is more accurate, it can also be used as 
autool.for désign and check of accuracy of analogs. .It is used here to 
examine the accuracy of A. It is also shown, as a by-product of this 
analysis that 
1) One of the assumptions upon which A is based is false, 
2) The BCG contribution of one group of the signals available from A 

has been ignored, 
This chapter also includes a new, simple, and more realistic model of 
groups of branching small vessels. The third part is, as mentioned above, 


a study of the contribution of the heart to the BCG. This problem is 
| 
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equivalent to an investigation of the motion of the heart during the 

heart cycle. Two such attempts have already been made: one by Dr. 
Noordergraaf et al (4), which is basically synthetic (and cannot, therefore, 
be easily and reliably used for routine simulation) - it is based on 
experimental data not available from A and on several simplifying 
assumptions. The other, proposed by Hooks (5), is analytical but has two 
serious drawbacks: it is oversimplified and completely omits two out of 
four essential factors which cause the motion of the heart. In the 
approach suggested here, physically more accurate assumptions are made and 
all four factors are considered. 

It is difficult to evaluate the accuracy of the approach proposed 
im the third chapter. Thére are three pose hbemaye of doing it. One is 
comparison with results obtained in previous studies (mainly (4) since’ (5) 
does not seem to be realistic). The correspondence is not complete; 

reasons for discrepancies are discussed in the appropriate place. Another 
iS a comparison of a real and simulated BCG - this is of limited value only, 
mainly because of the large variability of BCG's as mentioned earlier. It 
would be more valuable in some physiologically abnormal cases - but in 
these cases data necessary for calculations are not available. The third 
way is an analysis of the physical factors behind the simulated phenomenon. 
This too will be done later, 

There are two encouraging results: the contribution of the heart 
can be simulated under even more simplifying assumptions without affecting 
results noticably and some parameters, which are not known to an accuracy 
normally considered as sufficient, do not affect results substantially. 


The CGS system of units will be used throughout this thesis 
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because it is the standard system used in this field. 

To conclude this introduction let it be stressed here that this 
thesis presents methods rather than results. 

This work, although somewhat critical of previous research in 
this field, would not have been possible without data and results 
gathered and formulated by other people - mainly Dr. Noordergraaf and 
his students, 

It is the hope of this author that methods presented: here will 
eventually help in advancing the understanding of BCG and so lead to a 


broader application of this potentially very valuable method. 
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CHAPTER 1 


Baste equations and definitions 

In this chapter basic equations of BCG simulation are derived 
and used to divide the total BCG into several components which can be 
calculated separately. Although these equations are very simple, it 
seems that they have not been analyzed carefully enough by previous 
authors. It will be shown that this has lead to a conclusion which 
could possibly result in serious errors in simulation (see end of this 
chapter). 

To simplify notation only the head-foot (y) axis will be con- 
sidered. Similar equations hold in other directions. The reason why- 
this axis is discussed is that it is generally believed to be the most 
important one, in that y BCG's seem to contain more information than 
records taken in other directions. Also, it is easier, for technical 


reasons, to obtain reproducible records in this direction. 
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Notation (see also Fig. 1): 


Ay 


total mass of the subject. 

total mass of the frame of subject's body. 

mass of subject's blood. 

a fixed point in the external frame. 

a fixed point on the steady part of subject's frame. 

position of an infinitesimal volume of subject's body at 'rest' 
(corresponds to a fictitious state when the subject is not 
breathing and his heart is not moving but gives an output corres- 
ponding to the DC component of blood pressure) with respect to A. 
position of an infinitesimal volume with respect to the resting 
state at time t. (In general symbols without a prime have the same 
meaning as symbols with a prime but the reference point is B). 
position of the center of gravity (CG) of the frame at rest with 
BCSpecueLOmAs 

position of the CG of the, total volume of blood at rest with respect 
to A. 

part of the body volume corresponding to mm: 


part of the body volume corresponding to Mee 


Basic equations 


Essentially, BCG records reflect changes of the position of the CG 


of the subject's body due to breathing and circulatory events. In this 
6 y. 


thesis, only changes due to circulatory events will be analyzed, although 


BCG records are also influenced by mechanical coupling between the subject's 


body and the BCG bed, properties of measuring devices etc. These effects 
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can be represented with sufficient accuracy by a linear block in series 
with the source of BCG as described here. 

Let us examine the equation for the position of the CG of the 
subject lying on a freely’ floating BCG bed, assuming that at rest (as 
defined above) the CG is not moving. Since there are no external forces 
acting upon the body, the position of the CG w.r.t. A must be constant 
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here dm is an infinitesimal mass in the infinitesimal volume dV at rest, 
dm + dAm is the infinitesimal mass in this volume at time t. 
Introducing Me thewposition of A w.r.t. B gives (y' = y + Yan 


and Ay' = Ay) 


Yapl 4 (dm + dAm) + J (dm + dAm)] + i Ay.dm + i Ay.dAm 
MS Ve Viv, vitV_ 


re, [ y.dAm + [ y.dm = Const (1-2) 
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* Since the method of calculation of the center of gravity represented 
by equation (1-1) is not quite common a more detailed explanation is 


given in an appendix to this chapter, 
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+ mY, a MeYe = Const (1-3) 


Yap is the variable that we measure. By a proper choice of 
points A and B all constant terms can be cancelled out (they are of no 


importance anyway since it is the derivatives we usually measure). Let 


us assume that A and B have been chosen in this way. Then 
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S pail is due to changing volumes of blood in vessels and the heart 
- assumed not to be moving. 

yay is the contribution due to the motion of some vessels (e.g. the 
aorta) and their variable volume. 

ne is the contribution of the moving parts of the frame (e.g. the 


moving heart filled with variable volume of blood ete). 


It is not claimed here that these contributions are equally 
important - this remains to be evaluated, 

Let us now examine these rather obvious results from the point of 
view of simulation. It is generally believed that PER Tea constitutes 
by far the most Aer component of the total BCG. This is the main 
reason why attempts to simulate BCG by electrical analogs have been limited 
to simulation of this term (andin fact only a part of it: that caused oy 
the arterial parts of systemic and pulmonary circulation). Another reason 
is that this component can be simulated by linear passive elements without 
difficulties. One of the assumptions designers of this analog have made 
is that the choice of the reference point B is irrelevant. This is not 
correct when the quantity to be simulated is defined as it has been defined 
here (which seems to be the only reasonable way of doing it): assuming 
that pee represents the total BCG accurately enough means that all 
displacements are neglected. In particular the motion of the heart is 
neglected too. The heart then is assumed to be a pump which is not moving. 
But its volume is still changing. If this factor is not included in the 
model this leaves us only one possible choice of point B, namely B must 


be placed in the center of gravity of the pump. In the next chapter a 
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quantity representing the effect of a different choice of B is defined 
and evaluated. 

The ee ee and yee components defined above will be analyzed in 
detail in the following chapters of this thesis. Investigation of a fees, 
would require a more complicated model. It can, in fact, be analyzed once 


an acceptable way of analyzing vee has been found. 


APPENDUX £0..CH. 1 


The center of gravity of a body is usually calculated from one 


of the following equations: 


Vcekt)) sents) aloco oc). =< y(t)dm CLD) 
M M 
or 
eh 
Yog bt) peace yo (v,t)dv, (1-10) 
V 
where 


m = Jan = 0 (v)dv, 


y is distance to a fixed point on the external frame and v denotes 
volume. 
In (1-9) y is a funetion of time,and 6 constant, referring to 


a fixed element of mass dm = pdv. 





dial 


In (1-10) y is a constant (it is the position of the fixed 
infinitessimal valume dv) and p is a function of time. 

A third approach, which is a combination of the two described 
above, may be advantageous. It will be explained by an example, the 
case of elastic tubing in motion. 

Let us consider elastic tubing and a pulsating flow. Each 
infinitessimal Segment of the. tubing contains in general atime 
dependent volume of liquid and its distance to the reference point is 
a function of time. We can take these segments as the fixed elements of 


integration and write 


Yog 6b) -+ f yeeancey (1-11) 
V 


In equation (1-1) this approach is combined with the division of the 


considered system into two complimentary sub-systems denoted V, and 


b 
Vee In this case (1-11) becomes 
2 ie) Goi ORTGE (1212) 
Yoo = a by y(t)dm(t J y(t) dm 
Ms Me 


This is the meaning of equation (1-1). 
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CHAPTER 2 


Contribution of blood moving in the arterial part of the circulation. 


As stated in the preceding chapter, the component vetoed is 
believed to be the most important one. The venous part of the circulation 
is believed to have small effect because the pressure changes there are 
relatively small - and these changes are the decisive factor for BCG, as 
will be shown later. This chapter contains an outline of a new method of 
simulation of Sieod with two possible applications: 

Peeeitecan replace the analog*As-sirve. it fulfills the ee functions but 
more accurately and more economically. 

2) it can be used as a tool for the design of new analogs or a check of 
accuracy of existing ones. This last aneeaeia fe shown in this 
chapter and the accuracy of A is critically examined. 

This chapter also contains a simple and relatively realistic 
model of groups of small vessels, to be used in the method of BCG 
simulation described here. 

Only the y axis will be investigated for two reasons - no 
information necessary for the evaluation of contributions in other 
directions is available and the contribution of this part of the 
circulation (almost completely symmetrical in the head-foot axis) is 
believed to be negligible as far as the lateral (x) BCG (the only other 


BCG practically measureable) is concerned. 
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In the rest of this chapter the symbol BCG will denote only 


the contribution to be studied here, 
General considerations 


Several models, suitable for computer realization, will be 
described later in this chapter. All these programs are general in that 
the user can prescribe his own description of the system to be analyzed. 
All the models are based upon the assumption that the circulatory system 
(or, rather, those parts of it which are to be analyzed) is linear and 


that its distributed parameters are described by the following equations 


Ges 2) 
anes = Z)(s,2)F(s,z) (2-1) 
dF(s,z) PiCsez) 

dz rte Zits 2). re 


in the s-domain (Laplace transform w.r.t. time, zero initial conditions). 
Here P(s,z) is the Laplace transform of the pressure at distance z from 


the source, 


co 


P(s,z) = f e p(t, z)dt (2-3) 
0 
F(s,z) is the Laplace transferm of) the flew at z. 
Zoe are the distributed longitudinal and transversal impedances 
pet unit Length: 


An infinitesimal segment of vessel (parallel to axis y) of length 


dy gives the following contribution to the total BCG. 
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a(aeg(ed)=— F-[ly.ay, eae (254) 


(Integration with respect to t) 


This means that the BCG is a linear function of f. 
If we assume that the pressure (and, therefore, the flow) is 
harmonic 


t 


peehz) = P(w,z).et” (2-5) 


. (2-6) 


Ht 


£(t%z) F(w,z).e7- 


thae-means that the distributed BCG (CBCG') will have a similar 


form 


d(bee(t)) =, BCG! (a,2).e1 Oh dy (=D 


After substitution the following expression for the amplitude 
of the distributed BCG is obtained. 


: Pw, 2) 
BCGE(y, 2) == <y. 


Ss 


eis (2-8) 
JOZ 


Whatever the representation of this branching system the model has 
to be segmented and therefore the following equation has to be used to 
obtain the total BCG. 

BCGT(w) = y ABCG (2-9) 
all segments 
where ABCG is a BCG contribution of one segment at frequency w. 

This variable which is a function of frequency will be called BCG 

impedance, since it will be calculated from the assumption that the flow 
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Fw) =1 (2-10) 


Tf the actual flow is 


N N 
ACL). = >) Fonw).el"* . 2Reat() F(au) 23%) + FO) (2-11) 
-N “ n= 


(The Fourier expansion into a sufficient number of terms) 
the total BCG is 
N 


Begt(t) = 2.Real(\_ as 
n= 


BCGT (nw) .F(nw).e ) (2-12) 


a 


(The DC component is not considered since this constant is of no practical 
interest). 

The variable BCGT(w) thus gives all the information necessary to 
Obtain best for periodic input flows. This is quite sufficient since - 
beac lcitasare Of noypracticallinterest either. Also any complex of input 
flows can be approximated by Fourier coefficients as closely as desired. 


Tf p(t) is given instead of f(t) then 





Ao) 2 (2-13) 


and the impedance Z(w) is calculated by the program along with BCGT(w). 
To evaluate the effect of the choice of position of the reference 
point B upon the variable BCG'(w), let us shift B by a unit of length, 


This results in a change of BCG' by 


1 uy p.P(w) Ve 
Lhe) asta ey, dt Paulie 


and the effect on the total BCG is described by 
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QT (w) = [veer (2-15) 


“b 


Tf .B is shifted by L(cm) the new value of BCGT (BCGT) will.be 
BCGT(w) = BCGT(w) + L.QT(w) (2-16) 


If the correct value BCG(w) is to be obtained L.QT(w) must be subtracted 
from BCGT according to (2-16). 

The variable Q is evaluated by the program along with Z(w) 
and BCGT(w). 

The physical meaning of this variable can be understood from 
equations (2-14),(2-15). It corresponds to the surplus of blood 
accumulated in the part of the system being considered at the given 
frequency. When this variable is left out a change of position of the 
reference point B causes a change in the shape of the BCG frequency 
enaraeteristic., This, in turn, leads.to:a different shape of the BCG 
calculated from the frequency characteristic. But the apparatus is always 


measuring the same BCG. The variable QT(w)is, therefore, very important. 


Development of the approach 


The chronological development of various models will be described 
here - it shows which modifications (there are 8 different versions al- 
together) have been progressively introduced and why. 

The first decision that had to be made was whether the time domain 
or aie frequency domain should be used. It is evident from the previous 


section that it was decided to use the frequency domain. The reasons for 
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this choice are the following: 


1) 


2) 


3) 


4) 


the frequency characteristic of the system can be calculated once 

and for all and the desired studies can then be made via Fourier 
Meernictcnte: This saves a lot of computer time and even makes it 
possible to study BCG without a computer. On the other hand 
transients caused by changes of parameters cannot be studied since the 
source - the heart - is connected to the system by nonlinear elements. 
This OEE ETE is not important even when the method is compared 
with the analog A - the analog works 1000 times faster than reality 
and it is, therefore, hardly possible to study transient phenomena by 
its use (the system being a distributed one). The analog, however, 
fas the advantage that, to some extent at least, it generates its own 
input (flow). 

the frequency response is a good and accepted characteristic of the 
circulatory system. Several studies published so far make it possible 
to compare models and real systems via frequency characteristics (as 
far as the circulatory system itself is concerned). 

since the circulatory system is represented by a branching system of 
uniform tubes, each of which is an element with distributed parameters, 
the frequency domain is the only approach which allows a fast and 
general exact solution. It will be seen that this advantage over 
analogs is an important one. This feature of the proposed method also 
makes possible its use as a diagnostic tool. 

as a consequence of 1 and 2 a computer program using the frequency 
domain offers an efficient method for comparison of various physical 


models and their simplifications (provided they are linear) - also 
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possibly in the first phase of design of analogs. 


A. Original Model. 


This model’is called C when complex value of E (modulus of 
elasticity of vessel wall) is used, R when real E is used, 

This first version is a complete analog of A in digital form 
and frequency domain - i,e.\the topology and parameters are identical. 
Also, the system is divided into segments which are assumed to have 
lumped parameters. 

The system to be simulated is assumed to be described by 


distributed impedances 


1 oP 





OES re ary Coe 
' a oF 
2) (w) =- P(w)/s> C= 1e) 


z is the distance along the vessel measured in the direction 
from the source as mentioned above. In versions C and R, these impedances 
are lumped and a segment of length L is presented by overall longitudinal 


and transversal impedances (see Figs. 2 and 15). 


Z, = Ze (2-19) 


Zi det 2! [I Ce) 
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L = Inductance 
R = Resistance 
C = Capacitance 
ies 2 One segment of vessel as represented in (1). 


Expressions for 2) and zt are taken from (1) and are as 


follows d 
‘ a nip 
Au 0 
Zt Ab (2-21) 
Say Lite arid lk Glo] 
0 
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- eptsel wp Ar atieh 


a! (2-22) 
E ive area (x + says 
where 
aah 
Sree (2-25) 
here 


rt is the internal radius of the vessel when the internal pressure is 
equal to the DC component of the actual pressure, 

0 is the density of blood. 

w is the frequency w = 21.f 

h is the thickness of the wall of the vessel. 


E is its modulus of elasticity. 


Equation (22)contains a ratio of two infinite series (representing 
Bessel functions). A fast converging continued fraction expansion is 
given in (1). A few terms are then used to obtain a desired accuracy. 
The whole segment is then represented (in the electrical analog of which 
this version is a computer model) by the network in Fig. 2 of total 


impedance Z, and Zi 


x 

Although the accuracy of individual blocks is good, the large 
number of segments (over 100) and tolerances of the values of individual 
passive elements, as well as the effect of connections etc., limit the 
overall accuracy of the analog. The number of passive elements in 
individual segments was chosen such as to give 2% accuracy. In the 
computer program this was replaced by 12. 


Physical limitations of analog A are to some extent replaced by 


roundoff errors. These were minimized by the use of double precision 
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arithmetic for important variables. The accuracy of the computer program 
has been checked by a simple branching system resembling the circulatory 
system but having uniform parts. Exact solutions were then compared 

with solutions obtained by subdividing this system into more and more 
lumped segments. The convergence was excellent and thus it was concluded 
faat the program is accurate. Results obtained from this version of 
computer program were compared with published frequency characteristics of 
the analog A (see (1)) and the agreement of hemodynamic variables has been 
found good. This is shown in Fig. 6 and confirms the relative accuracy 

of the analog, at least as far as hemodynamic variables are concerned. 

No comparable characteristic of its BCG performance is available. The 
original version of the analog A was assuming real values of E, This 
assumption has then been replaced by complex F (physically elastic 

modulus of viscoelastic material of vessel walls). Both these alternatives 
were tested on the computer program and it was found that the second 
version of the analog gives results closer to those obtained from the 


digital computer. 


B. Distributed version of Original Model. 


This version is denoted by CD or RD for complex and real value of 
E Deon eetively: 

Unlike the analog model or the time domain approach on digital 
computers, the use of frequency domain makes it possible to obtain accurate 
results even for quite general linear systems with distributed parameters 
(by 'general' is meant generality of expressions for distributed impedances). 


It must be remembered, however, that even this solution will only be an 
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approximation because the continuous and general way in which physical ° 


parameters are changing must be approximated. 


When the basic equations 


dP’ 


de =- Z)+F (2-24) 
dF P 
7 oan (2-25) 


t 
are solved (e.g. (7)) (2) 2 are assumed to be constant throughout the 


individual segments), we obtain 


2 
Nea Cae ston Cyn lL) 
Be inCuUU eee (2296) 


in 2 ' 
Y Zap tgh(y-L) + YZ 


for s # 0 and 


ES L ' - 
Aes = Le Xe + LZ) (2-27) 
for s. = 10 
Here 
7 
Q 
; at (2-28) 
1B 
Tay 
Q P 
P@s) = Be: (cosh (7.2) + 7 aes sinh (y.2)) (2-29) 
out 
Here 


Ve a is the loading impedance at the end of the segment. 


Zs is the loading impedance presented by the segment and its load. 





Pe.) 


L is the length of the segment 
P(2£,s) is the Laplace transform of the pressure at distance % from the 


end of the segment. 


For the calculation of BCG the artery is assumed to have parameters 


indicated in Fig. 3. 
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In this figure 8 > 0 when Z Zn (measured from the source). 


< 
A 


For s = jw 


L iL. 
ABCG (w) sh BCG'(w) dg = ee [PG ,ary.ae (2230) 
0 A ars! e t 0 * 


but y = Me + 8.2 so that 





‘ hes 
ABCG(w) = cn con Cyt) bee —-— oS pe Ya) Cy se oe) ae 2—31) 
‘J out . 
Finally 
P eae 2 
ee eee ee el x 1 © 
ABCG (wu) Seer 7 tue + 6L)Y B.Z)}.sinh (yL) 
s is @ ike 
ee aS. 1 = 
+ mys ie BL)Z) BZ pe cosh (yL) + nena vot) (2-32) 


and AQ ( substituting 8 = O and vee ‘Bh 


mm; ue 388 (y¥L) + Z) {cosh ona) (2233) 


At this-point only was “it vealized that the choice of othe reterence 
point is important and the variable Q was introduced. The analog was then 
examined more closely and it was concluded that there are two other 
questionable points. Both concern the part of analog A representing the 
small vessels: 

1) In the analog small vessels are represented by resistors. This is 


adequate for hemodynamic variables (which was confirmed later when a 
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new model was tried) but perhaps not for BCG simulation. The reason 

is that resistors replace a system of elastic branching vessels by 

one element with zero time delay and no elasticity. The contribution 

of small vessels is (in the analog A) calculated from 
Bes 8b ie 

Oe eet) ty Pe) a (2-34) 

E.(2r + h)h 

with 

end £ = 16 x 10° (see (4)). 

V is the volume of the considered system of small vessels. 

The result ‘of this approach is that the contribution of small 

wessels tar from the reference point is suspiciously large. 
the analog does not take into account one important factor at all, . 
perhaps intentionally, but then the model lacks physical meaning. The 
model is at all terminations grounded by resistors. Some current 
(the electrical analog of flow), however, is flowing out or into the 
system through these resistors. Since the amount of blood in the 
system panne change - the physical system being closed - the model 
represents a system terminated by 'reservoirs', This fact is not 
considered at all. It was decided here to view the reservoirs as 
localized at ends of resistors representing small vessels (see Fig. 4). 
ties Pserathcer arbitrary, soutett iS cimple anda eat lease considers | 
the neglected effect. If f(t) denotes the pleodalon from the 


reservoir then its contribution to- the BCG is 


hPon (bys =e huCe(@) teu © & fe cery.at = fey.) .dt = ee eve 
Ss Ss 


(2-35) 
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Completion of the system by ‘reservoirs’. 
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ABCG(w) = ne 





(2-36) 


If the total effect of this factor were negligible, its omission 
would be justified, but this is not the case. 

For these two reasons it was decided that the model had to be 
modified. This is described in the next section. 


Cc. Model with modified representation of small vessels. 


This version is denoted by Cl and Rl for lumped version and CD1 and 
RD1 for distributed version. 

This new description is partly based on greatly simplified 
physiological data, partly on data used in the original model, and partly 
on the formula for distributed longitudinal impedance given earlier. 

Fig. 5 shows the geometric configuration assumed for branching small 
vessels. It is evident that, compared to the nonsymmetrical vessel anatomy, 
the model is very simplified. On the other hand the number of branches and 
groups as well as the ratio of volumes of successive groups were kept as 
close to reality as possible (for comparison see (6)). 

The symmetry of the system and some further assumptions make it 
easy to find parameters describing these systems of vessels. 


Let us assume that 
R =—— (2=37) 


Here RL is the resistance of one branch of the n-th group, N is 


the number of branches. Since the number of branches of the n-th group 
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is 1/N times the number of branches of the (n + 1)-st group, this means 
that the overall resistance of groups is constant from one group to 
another. 


Since 


(This can be obtained by substituting »’ = 0 into (2-21) section A).For 


one vessel, (2-37) gives 





g 2 
n 1 ntl 
—_—_—_ = 2- 
4 NA Ne, 
n ntl 
When we define 
Q Te 
ts n+1 a n+l Lis 
=) Ae? Cy or. ey, 
n n 
this means that 
4 
Cee INe (2-41) 
g te 


fe Rm) is the loading resistance of a branch of the (n-1)-st 


group (see Fig. 5) 


eM ote 420, (2-42) 
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since at the end R 


But according to (1) 


1 1 ik 
ee R —— R = ferue eae s (2244) 
Nace n2 2 wi Ry 


so that 


ot _ Suk M A 
R el en (2-45) 


ug 


For the first 5 sections (included to obtain a smoother transition from 


large to small diameters) let 


gf 
ae oe (2-46) 
i & 
n 
with 
Sa) =a a (2-47) 


here ty is ideally equal to the radius of the vessel immediately preceding 


the considered system of small vessels. Let 


C= Cae ts (2-48) 


then 


1 M 
LS = Be | Ge ap — Fat eet eee: 0) oo TR eee —na| (2-49) 
Oo 





sill 


By a similar reasoning we obtain the following Expression for the total 


volume of the sytem: 


DO 
: : 10 inp pee OD a! 
RL ovee (eeiteagcre. Vitlsauc)xels +8Na- & | (2-50) 
total 220 
1-WNe J 
i 
The ratio of total areas of cross sections of consecutive groups 
is 





nes 
age n+1 2 
shee tse a noe eae 
n r 
n 
Expressions (2-41)-(2-50) and data about rma ant and ey 


used in previous models are then used to obtain values for unknowns. There 
are 3 equations and 7 unknowns (5 2%sa,C),c M,N). It is therefore 
necessary to choose some of them arbitrarily - but as close to physical 


reality as possible (see (6)). It was decided to set 


Gee eC = 68 75 Cem 2D 
r Q s 

M r=. 6 

Nee LL 


A digital computer program was then used to find acceptable values 
for to by varying the value of a. 

The result of this approach is that the total volume and resistance 
remain almost identical (the difference is due only to the roundoff error 
introduced into the program) meh those used in the analog and the above 
mentioned versions. 


The 'reservoir effect" is calculated from outflows from the last 
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group according to (2-36) in the previous section. 

The model presented here is believed to be more realistic 
than those suggested by other authors. It satisfies the requirements 
for which it was designed as will be seen later in this chapter. 

It should es mentioned that a physical analog with these 
characteristics is practically unrealizable because of the enormous 
number of elements involved. The computer model, however, requires 
little computer time because of its symmetry. 


D. A simple model based on Poisseuille's Flow. 


This version is denoted LL for lumped and LD for distributed 
parameters. Only complex E was used for calculations. 
Welkovitz and Fich in (7) derived and used the following expressions 


for distributed impedances: 


zi = 22 + jw, 5 = Rete jusl (2-52) 
Tr Tr 

eee HE eT: z 

ciede 32a Giuc (2-53) 
jw2tr 


v in the first equation is apparently a printing error and should 
be replaced by u. These equations were, with some modifications, used by 
the authors to describe a part of the aorta. If these equations could be 
used instead of those used in previous sections, the result would mean a 
substantial saving of computer time or a significant reduction of the 
number of elements in the analog. For this reason, it was decided to 
compare results from other programs with two versions using these equations 


- one using lumped parameters, the other using distributed parameters. 
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In addition to the change made in the first equation mentioned above it 


was decided to replace C by 


a ere: + nya 


ial E h(2 r+ h) (2-54) 


This is the formula used in previous sections taken from (1). Its use 
makes it easier to compare results from various versions of the original 
program. The eee mentioned formula seems to be generally accepted as being 
better justified than the one used in (7). 

Only the modified model of small vessels described in the previous 


section was used. 


Results and Conclusions 


Frequency characteristics obtained from different versions are 
in Figs. 7, 8. It should be mentioned that 'force BCG' and its integrals 
rather than the more common acceleration BCG (acceleration of the BCG bed) 
and its integrals are used. This is because the author believes that 
dividing the force BCG by the mass of the tig (which is the transformation 
through which the acceleration BCG is obtained from the force BCG) 
introduces a somewhat modifiable scaling factor and thus obscures results 
of simulation. Only complex value of the modulus of elasticity E is 
considered. Real value of E was originally investigated by the dea 
mainly for the purpose of checking the accuracy of the analog. 

Examination of results leads to the following immediate conclusions: 
1) the overall hemodynamic impedance calculated by the program is quite 

close to that measured on the analog. This is true both for lumped 


(C) and distributed (CD) parameters versions. It means that from the 
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point of view of circulation itself the number of segments into which 
the simulated system is divided is sufficient. 

BCG frequency characteristics C and CD differ significantly in amplitude 
in the lower frequency range. This means that for BCG simulation the 
number of See should be increased. 

BCG frequency characteristics of the modified lumped (Cl) and 
distributed (CD1) model differ from C and D mostly at high frequencies. 
The modified version gives generally smoother results in agreement with 
reality. This difference is probably smaller when Cl and CD1 are 
compared with the analog A because the analog works 1000 times faster 
than reality and various unwanted effects can affect its frequency 
characteristic at higher frequencies. 

The relative contribution of small vessels is much smaller in Cl and 

CD1 than in C and CD. This can be seen in Fig. 9. This confirms the 
original impression that the contribution of small vessels is exaggerated 
in A. Evidently the way in which these vessels are simulated plays an 
important role. 

Comparison of normal BCG records obtained in the way described above 
(see Figs. 10, 11) and clinical records (e.g. (8)) shows that amplitudes 
obtained from versions with distributed parameters are closer to real 
values than those obtained for lumped parameters. The same observation 
is valid for results obtained from A (see (3)). This supports the 
earlier conclusion that for accurate BCG simulation the number of 
segments should be increased if the system is to be simulated by an 
analog. 


the simplified versions LL and LD give essentially the same results as 
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Cl and CD1 which means that a simpler representation of the hemodynamic 
phenomena is acceptable when details of local physical phenomena are 
not essential. 
It should be noted that the results presented do not cover the 
contribution of the pulmonary circulation which can be obtained by the same 


program with different data (available from (2)). 


Summary of Chapter 2 


A new approach to BCG simulation is presented. The author feels 
that it is more realistic, economical and easier to use than the existing 
analog model simulating the same phenomenon (i.e. a part of the human BCG). 
Results, once obtained, make it possible to simulate various BCG's even 
without the use of computers. This approach makes easy comparison between 
various linear physical models of the same physiological system and can be 
used as a tool when an analog of the circulation is designed or critically 
analyzed. 

The approach described in this chapter is used to investigate the 
accuracy of an analog built at the University of Pennsylvania and described 
dn) (152,3). sit issconcluded that although this analog’ simulates the human 
BCG quite well, it is not physically complete, simulation of small vessels 
is not quite satisfactory and the number of segments into which the 
circulatory system is divided is not large enough. A conclusion that this 
kind of simulation is not ideal is put forward. A new model of representation 
of small vessels both physically more realistic than the one used in the 


analog and easy to handle on a digital computer is proposed. Introduction 
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of this modification leads to a smaller contribution of small vessels 


and a smoother BCG in agreement with reality. 


“Appendix. to Chapter 2 


In the first place the organization of data giving information 
about the system to be analyzed will be described. Then the block diagram 
of the whole program will be presented and some blocks described in more 
detail. The program was written to be used via terminal under the MTS 


system implemented at the Computing Center at the University of Alberta. 


A. Data Description 


The system is assumed to be branching in the direction away from 
the source only. The following rules and definitions are accepted Sion the 
description of the system: 

1) 'Distal' nodes are nodes farther from the source. 'Proximal' nodes are 
closer to the source. 

2) All branches have both distal and proximal nodes (this is important for 
the terminal nodes). 

3) Every node and every branch has its number. 

4) Every branch has the same number as its distal node. 

5) Numbering starts from an arbitrary terminal node and proceeds in the 
way as one would proceed when calculating impedances of the system. 


These rules are illustrated by an example presented in Fig. 12. 
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Fig. 12 An Example Illustrating the Use of Some Rules from Section A. 


Block of segments (group of segments) number I is described 
by an integer vector B(I.J). 
B(I,1) = 0 for terminal nodes, 1 otherwise. 
B(1,2) = number of distal node of the I-th segment. 


B(L,3) = number of subdivisions into which the segment is-divided,. 


B(I,4) = 1 for the aorta, zero otherwise. (The reason for this distinction 
is that complex E is described by a different expression for the 
aorta than for other vessels). 

Parameters of the J-th subsegment of the I-th block are stored in 
the vector $(I,J,K). Components of this vector are the inner radius, wall 


thickness, modulus of elasticity, length and distance from the reference 


point of the segment etc. 


B. The Block Diagram of the Program and a Description of some of 


Its Blocks. 


The block diagram of the program is given in Figs. 13, 14. 
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1) Impedances of the lumped model are calculated from the configuration 


illustrated in Fig. 15. This leads to the following expressions: 


Lieve: Gi) ah He 


rae ~ Q We Teal, MORURE out (2-55) 
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Fig, 13 One segment of a vessel 
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2) A note on the order of calculations performed by the program. 
In each segment calculation proceeds in Sete et i out order’ 
At the beginning of calculations it is assumed that the outflow from 
the segment considered is F = 1. This means that in the case of 
lumped parameters Q and BCG are calculated from P and F at the distal 
end of the segment. For this reason calculations of “the contribution 
of the terminal group in the modified model of small vessels (which is 
very Ainoctant because of the number of branches involved) have to be 
performed using the assumption that parameters in this group are not 
lumped but distributed (unlike the rest of the system), since P = 
at the distal,end of this group, 

3) In part 3 of the program actual values of variables are found from 


assumed values (see previous section) from the condition of linearity: 


iP 
meals real 
= —— i 2— 
Per: iz P (2-57) 


Here ies is the flow through the considered node when the flow 


eal 


from the source is 
eee al (2-58) 


It is calculated in part 2 of the program. F is the flow through 
this node assumed in part 1. P is the pressure calculated in part l. 
A similar notation will be used for other variables throughout 


this section. 


(2-59) 
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BCG wal = LE Oe: (2-60) 


Q (2-61) 


real pres tae 
4) Modified model of small vessels. Using the same expressions as in 


section 3 and the fact that there are 11 branches and 6 groups 


P 
: 2g? aes (2262) 


reaL oon real, a P 


where Rene i ig the real outflow through one terminal node and F eal 


the total réal outflow from terminal nodes of the last group). 
In the calculation of BCG and Q the number of branches in each 
group is considered in part 2 of the program. For this reason 


BCG = F BCG aa(2=63) 


(2-64) 
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GHAPTER 3 


Contribution of the Heart 


This chapter examines the most significant contribution of the 
component Ne oor the contribution of the moving heart. 

It is generally accepted that the contribution of the heart is the 
second most important single component of the human BCG. Two attempts to 
estimate the motion of the heart have been made so far. One of them (4) 
has been made specifically for the purpose of BCG analysis,the other one 
for other purposes,but its possible application for BCG analysis has been 
expressed in (5). The first of these two studies is of a more synthetic 
nature - it is based on experimental data which are used to obtain the 
motion of the center of gravity of the heart via simplifying assumptions. 
The other is analytical in that it starts from a physical description of the 
system, makes some simplifying assumptions and arrives at the motion of the 
CG of the heart. Results of these two studies are quite different. Both 
these studies consider only the head-foot component of the BCG. the 
character of the approach presented in this chapter is analytical.) it 
uses some simplifications suggested in (5). It is, however, physically 
more complete. 

It is assumed in (5) that the heart is a sphere with a single 
outlet, the cross-section of which is equal to the total cross-section of 
vessels leaving or entering the heart (‘heart vessels' from now on). It 


is oriented in the head-foot (y) direction as is the direction of the blood- 
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flow. The bloodflow through this outlet is equal to the sum of bloodflows 
leaving or entering the heart. This means that the directions in which 
individual heart vessels are leaving or entering the heart are not 
considered, which is a very drastic simplification. As for other 
simplifications (the shape of the heart and its mechanical coupling to the 
chest cavity), these are made necessary by the lack of data and by physical 
complications which can hardly be overcome at this moment. The heart is 
assumed to be surrounded by a viscoelastic tissue which is, however, 
represented by lumped mechanical coupling (see Fig. 16). The tissue 
itself is enclosed in a cylindrical cavity (see Fig. 17). The only force 
(except for the force due to coupling and acceleration of the heart) 
considered is the force due to the flow through the only outlet. The 
dynamic equation governing the motion of the center of gravity of the . 


heart is then 


mY + dy +cy = Fe (3-1) 

with 
d 

vg Kin el a Ray 
here 
TM is the mass of the heart. 
c,d are viscoelastic constants (see Fig. 6)... (Their. derivations; ina) 

is described later. 

Fe is ie force due to outflows. 
Vet. is the relative velocity with which the blood leaves the heart. 


Comparison of results from (4) and (5) (Figs. 28,29) shows that they are 


completely different. 
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a. 


The first approach used by this author considered all six degrees 
of freedom (i.e. not only the displacement but also the rotation of the heart 
were simulated). It will be described here although it has been abandoned 
later for reasons that will be explained at an appropriate place in the 
Eext. 

Results obtained from this first model were used as a justification 
for simplifications introduced into models investigated later. 

The following description of the first model could, perhaps, be 
used as a starting point of another study. 

a) the heart is a spherical cavity the volume of which is given by the 
instantaneous volume of blood it contains and the volume of the heart 
muscle. It is realized that this assumption is rather arbitrary 
since the heart will be considered as a solid body coupled to the. 
external frame by mechanical springs and dashpots which neglects the 
motion of the vibrating tissue surrounding the heart. 

b) for simplicity it is assumed that the density of the heart muscle and 
of the blood is 1. 

c) the heart vessels considered are the same as the real ones. Their 
positioning on the surface of the sphere and directions in which blood 
flowing through them is leaving or entering the heart are chosen in 
such a way as to approximate the reality in an acceptable way. 

d) flows through various heart vessels are taken from (8) - presumably 
typical normal flows. 

5 the coupling of the heart to Fe walls of the thoracic cavity (which 
is assumed to be at rest - since the frequency of breathing is normally 


much lower than that of the heart, this seems to be an assumption which 
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should not affect results significantly) is viscoelastic. (In (5) 
only forces are considered). 

(£) coupling parameters are calculated from assumed viscoelastic parameters 
of the tissue in the same way as in (5) for forces and from simple 


assumptions for moments. 


It is obvious that these assumptions simplify the reality 
significantly. The main reasons why they were accepted are, Ail. Tene 
complexity and poor documentation of the real system. Without drastic 


simplifications simulation of this problem would not be possible at all. 


a) Parameters 


The mechanical coupling of the heart to the chest wall is assumed 


fo be aS in Fig. 16. 


Chest Wall 





Fig. 16 Mechanical coupling Heart - Chest Wall 


Index I = F is used for force coupling, M is used for rotational coupling. 
The geometry assumed in (5) is accepted as the basis: for the cal- 
culation of constants ¢ and d (Fig. 17): Constants for force coupling 


were derived from measurements published in (9)and the assumption that the 
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Fig. 17. Geometry assumed in (5) 





54 


heart is supported by a cylinder of a homogeneous tissue, the cross-section 
of which is equal to the cross-section of the heart. On the basis of these 
assumptions two sets of constants are calculated in (5) - one for head-foot 
and one for lateral coupling. These assumptions and constants were 
accepted here aie abe only change that constants were assumed to be the 
same in all directions. The reason for this simplification is that the 
geometry is so simplified that making the distinction between the two 
directions is hardly justified and also because the way in which coupling 
constants are calculated for momentum makes a similar assumption. For 


these reasons instead of 


5 sf 4 
Coy = 6.346 x LO doy =e 12929 ex 10 
Spire 627 SaeLOR Wi See ay 
iia eat Ps FL ; 
(head-foot and lateral) 
= by Sar Baethae ee 
Sigil Sk Dewees 
fi 
eh i i eee 


was assumed, 
To derive Cape dis the geometry illustrated in Fig. 18 was accepted. 
It is assumed that the heart is a cylinder coupled to a fixed 
outer ring by a viscoelastic material. (In the dynamic equation of the 
heart the mass of the surrounding tissue will again be neglected). Let 
us assume that under small deformation planes perpendicular to the axis of 


rotation remain planar and perpendicular to the axis of rotation and that 






d lianeire a 
Ye 


om of tandw Jos 









ows ozone 

i aot web 

vod be 

he .woleneath Lip ae 
— 

‘¢ notte : 
aa 


7 H Jie elite 2 


334 % 


} | 
' 
; 
~ 
igi le 
= 
{, 
} 


' mMeevaon ter LiL. Fa) As YG ‘yg with of 
+ 02) beabyaes holly: 
= 


e @t° aonaifredy, 1663 bootie tae 2 
n 7 - 


CE Ne Aid al rant 


| 
in? a 





55 





Chest Wal CZ. 
“ > 
RHdR-— 5 | 
o(R) + do 
| | AVR 6) 
S oe es 
, Ry 
M 
M 


Fac. 18 Derivation of constants C,; qe 


radial lines in these planes remain radial under the considered deformation. 


Then (see Fig. 18) 


a " M 
Tooth ee MCR) = Mo Ft eauy (3-3) 
dg = 5 with ¢ = ——- (3-4) 
2 la oe 
R R 
Ge oy ou M ») 
ies) a= | Ae Re TS, oy 
R 
a u 
here 
G is the modulus of rigidity of the tissue. 
I. ts Poisson's ratio which for incompressible material 1s .a (this is 7a 


rather arbitrary assumption but the effect of change of parameters 


proved to be small) 
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€ is shearing stress, 


At the same time we have (see Fig. 16 for coupling) 





or 
M 
eu r (sai) 
so that 
M 2 TRE 
Pesci = oe) (3-8) 
om M eee. ee ee 
e272G R Ry 
1 
here R) is ethe fradius of pthe heart, R, the radius of the chest cavity. 
3) : 6 
Roreeli=<L0ry 2-— 110; Ry = 5, R, = 12 we obtain Cy Ie? Tel Ome 
From 
2 
2mR)¥ 
che) CA gers eee (3-9) 
F a 1p 
Z je TAL 


(here Ro is the radius of the chest cavity). 


(This equation follows from the assumption that the heart with a 
total cross-section mR is supported by a cylinder of tissue (viscosity v) 
i 
' eee 
of effective length 9 5 Rj). 
and 
da = po 22 10° 
(calculated in (5) from data published in (9)). 
we get - 


v = 800 


and 
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As stated before the same constants are assumed in all directions. 


Summary of parameters: 
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b) Dynamic Equations 


The heart is assumed to have 6 degrees of freedom. The motion of 
the heart will be fully described by three equations for the displacement 
of the center of gravity of the heart, and three equations for rotation. 
These have to contain or be accompahied by equations for forces and moments. 
Let us start with equations governing the motion of the center Of eravity, 


of a body with changing mass. It is well known that the following equation 


dm. 


_ ak _ 


Spe u outlets “rel, dt eet (3-11) 


describes the displacement of the CG of such a body. 


r is the displacement of the center of gravity of the body. 


ms is the time dependent mass of the body. 


rely is the relative velocity with which the mass is leaving the body 
through the I-th outlet with the rate of mass flow. 
Fe is the vector sum of external forces acting upon the body. 


A similar equation can be derived for rotation: 


Let L denote the angular momentum of the body: 
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La= ee nt RY dm = voliees Te x dm (3-12) 


have 


Les Jia twith oJ = orc (3-13) 


w is the angular velocity of the body. 

The equation for the change of angular momentum during an 
infinitesimal time interval dt under the action of external moment ue 
is 


L(t + dt) + L, (t Je flisgh L(t) = Moat (3-14) 


(subscript H for the hearr, B for the blood that left the heart in the 


amnterval dt) or in our case 


(J+ 4) Ww + dw) + oe x ie Sed M @G=15) 


here coe is the position vector of an outlet with respect to the center 
of the sphere, u the velocity of blood leaving the heart. We will assume 


only one outlet for simplicity. This leads to 


dJ - du, - = dm) 
dt” ers dt a rout ar dts ee aly 
or 
Se a oe Sac Thee 
Jw + w - + eee. Oleg Ts Mest (3-17) 


when we replace u by 
(3-18) 
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38 rout : 3o51 dt rs (3-19) 


External Forces 


We have two kinds of external forces - one is due to the mechanical 
coupling, the other is due to the fact that the heart is connected to 
heart vessels which pull at the heart because of the forces generated in 
them by the Shayne direction of/blood flowing in them. This concerns 
mainly the aorta and the pulmonary artery. 


Forces due to coupling are described by 


LE Tee shes (3-20) 


A similar equation is used for moments, 

The last component of the external force to be considered here is 
the force developed by the changing direction of flow in the heart vessels 
which is transmitted to the heart through the walls of these vessels. 

Only the force developed in the aorta and in the pulmonary artery will be 


used, the others being considered negligible. 


Aw thes Aorta 

It will be assumed that the shape of the aorta is circular right 
from the aortic root. Since a part of the blood flowing through the aorta 
is leaving it approximately at the highest point of the aortic arch (see 
Fig. 19) an outlet will be assumed at this point and the amount leaving 


here will be taken as 


ae Oey, Seott) (3-21) 
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Aortele Arch 


Fig. 19 Tire Aortic Arch 


here ae is the bloodflow at point 2. 
f is the bloodflow at point 1 immediately before the outlet. 
Evers taken top beeequal ito -.4. 
Te is vealizedwthat this is’ 4 simplification. 
An infinitessimal mass dm moving along an arch of radius R with 


velocity v causes a force 


2 
a = (E22) 


which is pulling radially from the center of the arch. 
Let us assume that the velocity of flow at time t and point (R,$) 
(polar coordinates) (see Fig. 20) is (assuming an aorta with uniform 


parameters) 


vit, oa v(t = Rica) SOB AG RGpaa)y for ¢ < 1 (3-23) 
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Hen ee heen Rr) 1 IS RE: Os (3-24) 


where ons the velocity of transmission of a pressure pulse. Let us 


take ae 500 cm/sec. 






Aortic Arch 


Fig. 20 The Aortic Arch 


We have 





= SB oncepapenee ME RK, (3-25) 


where A is the area of the cross-section of the vessel, assumed to be 
constant. 

For the two components of the total force in the plane of the 
aortic arch we have (assuming p = 1) 


B . B 
F = [ dF. = He eR OK sind do (3-26) 
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liane. 
Been Pdhie sion eveCe yA tes} do (3-27) 
J Qa 7 J Qa 
Expressions (23,24) are substituted into (26 rand (27 )eandrthe 
resulting integrals evaluated numerically. The force Fa is then 
decomposed into its components Fos F wecordimeito, Fig. 21, (Z is the 
vertical axis, X horizontal oriented to the right hand side of the 


subject). 





Fig. 21 Decomposition of F 


Heart H’ 


The next step is to obtain the force acting upon the heart. To 
do this it is assumed that the physical situation can be simplified to 


that illustrated in Fig. 22 (zero moment) 





Fig. 22 Decomposition of pat 
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The component Fy is then obtained from 
(3-28) 


Fy Fo are assumed to act directly upon the heart. 
Bee tbe Pulmonary Artery 
With the notation of Fig. 23 we can derive the force due to the 


change of direction of flow in the following way: 
dm.v +E ..dt = dm.v_ (t + dt) = 0 (3-29) 
ai yP Me 


dt is the time in which the velocity of dm changes from ‘ CORZeron wile. 


changes direction by 90 degrees. This gives 


dm 


F Ve B30 
yp y dt ( ) 


This force acts upon the heart directly and is parallel’ to, the 
WV axis. 


We assume again 


Q 
RAG meg Ay) aes WAGE he 29). (3-31) 
Pp 
and similarly 
cL ae 
aoere he - ~,0 : (3-32) 
D ; 


2 is the length of the main pulmonary artery. 
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dm at time t 


Pulmonary Artery 


Heart 


Fig. 23 Flow in the Pulmonary Artery. 


This concludes the description of forces acting upon the heart. 
Let us note that because of the effect of smallness of distances with 
respect to c the effect of forces due to outflow and forces due comune 
change of direction of the flowing blood in the y direction practically 


cancel one another except for the first part of ejection - see Fig. Zire 


~C. Physiological and Anatomical Data 


It was not easy to gather the necessary data and some new assump- 
tions had to be made. Flows in various heart vessels (see Fig. 25) are 
# 
from (8,13) and represent a normal heart. Data about the mean volume of 
atria and ventricles as well as the mass of the heart muscle were taken 


from (15). The following table (Table 1) summarizes cross-sections 


and positions of outlets and directions of flows, which are deduced from 


a 


* Flows measured in dogs are extrapolated to man 
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PA Pulmonary artery 
A ‘Aorta 


VC, Vena Cava Inferior 


VC, Vena Cava Superior 
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illustrations in various anatomy textbooks. Positions of outlets are 
given as positions of points on a sphere, the diameter of which is the 
mean diameter of the heart corresponding to its mean volume. Directions 
of flows are given by end points of vectors of unit length (this is 


convenient for calculations). 


(axis z is perpendicular to the sagital plane, y is head-food). 


Table £. Some anatomical data. 


vena cava pulmonary 
' veins 


SUDEEIO“N SMnLerior 


Position of 
outlets 


Direction 
of flow 





To solve the described system of equations a computer program 
using the HPCG subroutine from the scientific subroutine package of IBM (14) 
was written. In order to include the effect of rotation of the heart, 
the motion of three points (at time 4 = 0 unit vectors in the x; vy andez 
directions) connected to the moving heart is calculated. We thus obtain 
a set of 21 simultaneous differential equations: 6 for the position and 


velocity of the CG of the heart, 6 for angular position and velocity 
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describing the rotation of the heart and 9 for the position efminree 
reference points. The last three vectors were continuously normalized 
in order to prevent the effect of numerical error from distorting the 
mutual positions of these orthonormal vectors. The details of this 
program will not ne described. 

When this set of equations was solved it was clear that the 
resulting rotation (as described by the three reference vectors) is 
relatively so small (the largest component — originally.ot Length 2, 
has never been changed by more than FOS)imethats in. view ofvother 
simplifications it seems admissible to consider only the motion of the CG 
of the heart instead of the whole rotating sphere. This is very important 
since the number of equations to be solved simultaneously is reduced 
from 21 to 6. Since these equations are decoupled (the effect of variables 
from one component of the three vector components does not influence 
equations for the two remaining components) we can solve just 2 or 4 of 
these equations if we are not interested in all components (which is 
usually the case). 

Only displacement (in the y axis) was investigated from this point 
on. As it was decided to simulate the effect of several factors (to be 
discussed later) separately it would have been expensive to do it on the 
complete model and the simplified model neglecting rotation was simulated. 
For this reason only results obtained from this simplified model will be 
shown and discussed here. 

Another change introduced at this point concerns the value of 


mechanical constants. The way in which constants, c,d, were obtained in 
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(5) is not very reliable because of the uncertainty of the value of the 
complex modulus of elasticity. (10) contains values of the natural 
frequency and damping of the heart inside the thoracie=cavity (f = 8 cps, 


6 = .25). If it is accepted that the system can be described by 


mY + dy + cy = Baw (3-33) 
then free oscillations are given by 
at 
y(t) =ae cos(wt + >) (3-34) 
with 
d 
g =- > (3-35) 
2m, 
and 
oe) 4mie - a ee (3-36) 
2m, V H r 
Since 
y 
6 = 2 (3-37) 
71 


(yy and y, are amplitudes of two waves separated by time = 
m > f. S$ can be used to calculate c and d. The values obtained 


and used from here on are 


(a esd (8) ne 10° iG heh a ee 10° 


When assumptions described above (except for rotation) were used 
to simulate the y BCG, results in Fig. 27 were obtained. It was then 


decided to try to incorporate the effect of the changing distribution of 
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blood in the heart and the shape of the heart itself during the heart 
cycle. It seemed possible to do this by changing equations for the 
coupling forces since the tissue is touching the. surface of the heart, 
assumed to be spherical, the force generated by the mechanical coupling 
is related to the Re keion and velocity of the center of the sphere. 
This is a good opportunity to incorporate the changing distribution of 
mass inside the heart. 

Let us of ee the following model of the heart: 

The spherical heart is hollow and divided into compartments - let 
us consider only two compartments for simplicity (one has a volume equal 
to the combined volume of both atria, the other has volume equal to the 
combined volume of the two ventricles). Let us assume that these volumes 
are constant (equal to the mean volumes) but that the mass contained in 
them is changing according to real volumes of the respective chambers, 
From this model the instantaneous position of the heart's CG with respect 
to the center of the sphere can be easily calculated. If 2 denotes the 
position of the CG of the heart along its longitudinal axis w.r.t. the 
center of the heart, Lis Le the fixed distances of CG's of the two 
compartments to the center and mm), ™ masses of blood in these compartments; 


Zz 


mM pthe total mass of blood inside the heart, we have 


me ifm. UY amet (3-38) 


Since the two compartments have nearly equal mean volumes (their 


centers of gravity are assumed to lie on the axis of the heart) we have 


¢ ME 1R =e (3-39) 
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SOmenae 


m2 = (m, - m5) &q (3-40) 


The force due to mechanical coupling is then calculated from the 


position of the center of the sphere as 
Bepeec(e —eh)jee a(t. 2) (3-41) 


To recapitulate: the equation of motion contains coordinates of 
the center of gravity of the heart, which is not identical with the 
geometric center of the spherical heart. The force due to coupling is 
calculated from the position of the center of the epneres which is obtained 
from the position of the CG calculated from the equation of motion saul Soar 
the calculated position of the CG w.r.t. the center of the heart. 

No experimental values have been found for the time-dependence of 
the flow between the atria and the ventricles. For this reason a calculated 
shape of the time dependence of the volume of the ventricles has been 
taken from (11) and the volume of the atria calculated from this and flows 
through the heart vessels (Fig. Os 

Results obtained from this modified model are shown in Fig. 24% 
which also includes the acceleration due to flows and the combination of 
these two factors. 

The following modification was then considered: as stated above, 
the heart does not really behave as an isolated single mass but rather as 


a mass surrounded by a continuum moving with it. It would not be 
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A Noordergraaf (4) 
B CG Contribution 
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practically possible to investigate this situation in this form. It 
seems, however, that the model could be made more realistic by lumping 


the continuum in a way illustrated in Fig. 31. 





Fig. 31 A Symmetrical Model of the Surrounding Tissue. 


Pnyorcablysthis means) enatethne thoracic cavity (or that partlor Lt 
which is assumed to be affected by the motion of the heart) is divided 
into sections, which are coupled to one another by springs and dashpots. 
The coupling was chosen as symmetrical and constants such that the 


overall spring and dashpot constants remained unchanged (Fig. 32). 


me 





Fig, 22 Equivalent Springs. 
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Cay Ca (3-42) 
d,=d (3-43) 


Results obtained from this model are in Figs: 27, 50% 
It was decided to try several sections to see the effect of 
lumping. For this purpose the situation illustrated in Fig. 33 was 


simulated, 





Fig. 33  Lumped Continuum, 


In this case 


o=5e (3-44) 
dues ad (3-45) 


The last modification was that the force due to the changing 
direction of bloodflow was taken as acting upon the section immediately 
above the heart (Fig. 34) and not the heart itself. Results are jin 


Fig. 27 which also includes the contribution of the heart due to outflow 





ree 


only. 





Fig. 34 Lumped Continuum. 


Discussion of Results 


Before results obtained from the described models are compared 
and discussed, some general notes must be made: 

Probably the only part of the BCG contribution of the heart that 
can be compared with experimental results is the H wave. This is because 
the H wave occurs around the beginning of the ventricular systole (12) 
before any significant effect of the remaining part of the circulation 
can be observed and also because the contribution of the normal heart 
in other parts of the BCG is relatively small. When the result published 
in (4) (Fig. 29) is examined from this point of view and compared with 
normal BCG's (11), it must be concluded that the amplitude of the H wave 
obtained in (4) is about twice the normal one (taking into account the 


variability of the H wave in normal BCG's). 
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Let us now compare results obtained from the model, considering 
forces due to outflow and reactions but not the effect of the CG of the 
heart, with results from (4). It can be seen from Fig. 28 that the 
timing ae the H wave and its orientation agree. The rest of the two 
results do not agree quite so closely but, as already mentioned, it is 
difficult to decide about the contribution of the heart except for the 
H wave. The amplitude of the H wave as calculated here is smaller than 
that obtained in (4) and also smaller than amplitudes normally observed 
(between 30,000 and 50,000 dynes). Reasons for the last discrepancy 
could be in the shape of flows used for calculations (the rate of change 
of flow plays an important role as can be judged from Fig. 24), in the 
accepted geometry of the aortic arch (larger distances lead to longer 
delays between the two forces in Fig. 24) (this seems to be in agreement 
with the fact that people with otherwise similar BCG records may have 
quite different amplitudes of the H wave), in the assumed dependence of 
flow velocity upon the distance from ate heart, the way in which the 
reaction from the aorta and the pulmonary aorta are transmitted to the heart 
and the way in which the continuum was represented by segments in our model. 
Also the effect of the eee Seen has notebeen includedsyet...eLE 
should be noted here that the contribution of the heart calculated in (4) 
has been redrawn for this thesis from a small picture with a different 
time scale so that some differences can have origin in inaccuracies 
resulting from this approximation. 

Results obtained from sa annniGne made about the effect of the 
changing blood distribution inside the heart cavity are not acceptable 


when the criterion of comparison with the H wave is used. 
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This may have two reasons: 

1) assumptions used for the calculation of the CG of the heart are not 
realistic (for example the real distribution of the heart muscle 
is not considered at all). 

2) knowledge of time courses of ventricular and atrial volumes used for 
the calculation of the CG of the heart is not good (this is an important 
factor because we are differentiating these variables twice and thus 
local details are significant). 

The assumption of ‘lumped continuum' not only seems to be realistic 
but also improves results of calculations. 

Results from (5) are summarized in Fig. 28. They will not be 
discussed here since the author thinks that they were obtained from an 


overly simplified description of the system. (see page DO) 


Summary of Chapter 3 


This chapter contains the first complete enumeration and discussion 
of factors associated with the motion of the heart - forces causing it and 
reactions. It is suggested that a model considering only displacements 
and neglecting the rotation of the heart is sufficient for the simulation 
of the BCG contribution of the heart since the effect of rotation is 
small. 

Factors affecting the motion of the heart are: 

1) forces due to the blood leaving the heart and entering it. 
2) forces due to mechanical coupling: heart - surrounding tissue - chest 
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3) forces due to the changing directions of flow through the eee 
arch and the pulmonary artery. 
4) forces due to the changing shape (and position of the CG) of the 
heart. 
S)artorces,due Poetic mechanical coupling with the tissue sur rounding 
the heart and moving with it as a continuum. 
Equations describing factors 1 to 5 are derived and used for the 
simulation of the motion of the heart. For the simulation of factor 4, 
a simple model making use of the coupling of the surface of the heart 
with the surrounding tissue is suggested. Factor 5 is simulated by a 
lumped representation of the motion of the continuum. 
Individual factors and their effect on the motion of the heart are 
examined separately (factors l and 3 are considered together) because 
the differential equations representing the motion of the heart are linear. 
This makes it possible to estimate the validity of description of 
individual actors. 
It is concluded from results obtained from this simplified approach 
that all factors, except for factor 4, are described adequately and that 
a reformulation of their representation will not be necessary in future 
research. Factor 4, on the other hand, is not described adequately and 
more research will be needed to simulate its effect. This author believes, 
however, that the method of ‘representation of this factor as proposed here 
is valid and that only a more sophisticated geometrical description and 
a better knowledge Eerene hemodynamic variables affecting 1ts influence 
on the motion of the heart will be needed to obtain realistic results. 


This will require a complete model of the circulation. For this purpose 
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analogs seem to be better suitéd than digital computers. The amount 

of experimentation involved in simulation could be too expensive if a 
dtoital computer were used. In) fact, once 4 linearized description 

of the problem at hand is accepted it is no longer justifiable to use a 


digital, computer for its)simulation. 
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Summary 


In the first part of this thesis a new method for -the 
calculation of the largest part of the human BCG (due to the arterial 
part of the circulation) is proposed. It is easy and economical to use. 
It can also be used as a diagnostic tool for the design of analogs or 
for checking the accuracy of existing analogs. It is used for this 
purpose here and an analog of the human circulation is evaluated from 
the point of view of BCG simulation. A new way of representation of small 
vessels for BCG simulation is suggested. 

In the second part (chapter 3) the contribution of the heart to 
the total BCG is simulated. Two methods for the solution of this problem 
have been published - (4,5), one of them (4) cannot be used for direct 
simulation (by an analog or digital computer). The other is over- 
simplified and incomplete. The approach suggested here is complete in 
that it considers all physical factors causing. the motion of the heart. 
More research will be needed to describe more accurately those of the 
considered factors which are not described adequately. It has been 
concluded that a relatively simple approach can be used to simulate the 
contribution of the heart mohtonuby neglecting the rotation of the heart. 
This is very important since nonlinear simulation is difficult by analogs 
and a detailed study of various factors (e.g. the calculation and effect 
of changing ventricular and atrial volumes) impossible unless a model of 
the complete circulation is used - and this would be very expensive on 
digital computers, 


It is the impression of the author that all components of the 
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human BCG have now been adequately described and that the simulation 
of the complete BCG is now possible. This is very important since 
detailed analysis of BCG records is the first condition for its 


complete understanding. 





(1) 


(2) 


(3) 


(4) 


(5) 


(6) 


(7) 
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